install.packages("foreign")
install.packages("RSiena")
install.packages("sna")
install.packages("network")
install.packages("statnet")

library(foreign)
library(RSiena)
library(sna)
library(network)
library(statnet)

#note: network data are symmetric
#basic scripts were based on https://www.stats.ox.ac.uk/~snijders/siena/RscriptSienaSymmetric.R 

#set your own work directory
setwd("~/Desktop/Research Wang/II UIA network/data")

wave1<-read.csv(file="wave1.csv",header=F,stringsAsFactors=T) #raw network data, in adj matrix format
wave2<-read.csv(file="wave2.csv",header=F,stringsAsFactors=T) 
wave3<-read.csv(file="wave3.csv",header=F,stringsAsFactors=T)

#read network files
wave1<-as.network(wave1,directed=FALSE, loops=FALSE, bipartite=FALSE, matrix.type="adjacency") #year 2002
wave2<-as.network(wave2,directed=FALSE, loops=FALSE, bipartite=FALSE,matrix.type="adjacency") #year 2007          
wave3<-as.network(wave3,directed=FALSE, loops=FALSE, bipartite=FALSE,matrix.type="adjacency") #year 2012  

wave1 #n of ties 2949
wave2 #n of ties 3230
wave3 #n of ties 3375

#calculate network attributes for three waves
network.density(wave1) #0.1528
network.density(wave2) #0.1673
network.density(wave3) #0.1748

#transitivity
gtrans(wave1, mode="graph", measure="weak") #0.3322796
gtrans(wave2, mode="graph", measure="weak") #0.3633601
gtrans(wave3, mode="graph", measure="weak") #0.3925971

#centralization
centralization(wave1, FUN=degree, mode="graph",cmode="degree") #0.7992412
centralization(wave2, FUN=degree, mode="graph",cmode="degree") #0.7948456
centralization(wave3, FUN=degree, mode="graph",cmode="degree") #0.7717949

#return all network files to matrix 
net1<-as.matrix(wave1)
net2<-as.matrix(wave2)
net3<-as.matrix(wave3)

#transform the network matrix data into Siena objects 
uia.t123 <- array(c(net1, net2, net3), dim=c(197, 197, 3))

# Create a Siena network object with sienaNet()
uia.net <- sienaNet(uia.t123, type="oneMode")

#read node attributes 
gdp <- as.matrix(read.csv(file="gdplogged.csv", header=T, stringsAsFactors=FALSE)) 
trade <- as.matrix(read.csv(file="tradelogged.csv", header=T, stringsAsFactors=FALSE)) 
mobile<- as.matrix(read.csv(file="mobile.csv", header=T, stringsAsFactors=FALSE))
internet<- as.matrix(read.csv(file="internet.csv", header=T, stringsAsFactors=FALSE))
pf<- as.matrix(read.csv(file="pf.csv", header=T, stringsAsFactors=FALSE)) #press freedom
hr<- as.matrix(read.csv(file="hr.csv", header=T, stringsAsFactors=FALSE))

#create country variables as changing co-variates
gdp.vc <- varCovar(gdp)
trade.vc<-varCovar(trade)
mobile.vc<-varCovar(mobile)
#class(mobile) <- "numeric" 
hr.vc<-varCovar(hr)
#class(hr) <- "numeric" 
#class(internet) <- "numeric" 
internet.vc<-varCovar(internet)
class(pf) <- "numeric" 
pf.vc <- varCovar(pf)

#dat.1 <- sienaDataCreate(uia.net)
#eff.1 <- getEffects(dat.1)
#effectsDocumentation(eff.1)

#jaccard coefficient wave1-wave2 .649; wave2-wave3 .679
#good for siena

#create the data for siena modeling
dat.1 <- sienaDataCreate(uia.net,gdp.vc,trade.vc,hr.vc,mobile.vc,internet.vc,pf.vc)
eff.1 <- getEffects(dat.1)
effectsDocumentation(eff.1)

#install.packages("xtable")
#library(xtable)

#set effects
eff.1 <- includeEffects(eff.1,inPop) #h1, degree of alter
eff.1 <-  includeEffects(eff.1,transTriads) #h2
eff.1 <- includeEffects(eff.1,altX,interaction1="gdp.vc") #h3
eff.1 <- includeEffects(eff.1,altX,interaction1="trade.vc") #h4
eff.1 <- includeEffects(eff.1,altX,interaction1="pf.vc") #h5
eff.1 <- includeEffects(eff.1,altX,interaction1="internet.vc") #h6
eff.1 <- includeEffects(eff.1,altX,interaction1="mobile.vc") #h7
eff.1 <- includeEffects(eff.1,altX,interaction1="hr.vc") #h8
eff.1

#try one model type = 2, indicating that the network is symmetric
mod.1 <- sienaModelCreate(projname='uia_model1',modelType=c(uia.net=2))
result.1 <- siena07(mod.1, data=dat.1, effects=eff.1)
result.1

#note that the results were generated from simulations, and thus the parameter coefficients may be different every time the model is run. 
#However, the coefficients should be similar in values.


